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ABSTRACT 

We present a new approach for simulating the atmospheric dynamics of the 
close-in giant planet HD209458b that allows for the decoupling of radiative and 
thermal energies, direct stellar heating of the interior, and the solution of the full 
3D Navier Stokes equations. Simulations reveal two distinct temperature inver- 
sions (increasing temperature with decreasing pressure) at the sub-stellar point 
due to the combined effects of opacity and dynamical flow structure and exhibit 
instabilities leading to changing velocities and temperatures on the nightside for 
a range of viscosities. Imposed on the quasi-static background, temperature vari- 
ations of up to 15% are seen near the terminators and the location of the coldest 
spot is seen to vary by more than 20°, occasionally appearing west of the anti-solar 
point. Our new approach introduces four major improvements to our previous 
methods including simultaneously solving both the thermal energy and radiative 
equations in both the optical and infrared, incorporating updated opacities, in- 
cluding a more accurate treatment of stellar energy deposition that incorporates 
the opacity relevant for higher energy stellar photons, and the addition of explicit 
turbulent viscosity. 

Subject headings: hot-Jupiters, atmospheric dynamics, radiative transfer 



1. Introduction 



Close in gas-giant planets have by now become a familiar part of the growing family of 
extrasolar planets. Their short period orbits and proclivity for transiting has made them 



the target of numerous observational campaigns and our knowledge of their structure and 
composition has increased dramatically over the past few years. The intense irradiation 
they receive from their host star dominates the energy budget of the atmosphere and drives 
supersonic flows unlike any seen in our own solar system. However, despite their prevalence 
and important role in constraining a wide range of planetary models, fundamental ques- 
tions about the dynamical behavior of their atmospheres remain, crucial for interpreting 
observations. 

There have been a number of different groups working to understand the dynamical 
redistribution of energy in close-in irradiated planets utilizing a wide range of approaches. 
Because of the breadth of these methods, interpreting the results and understanding their 
implications can be quite confusing. Det ailed comparisons of many of the models in the 



literature have been presented previously (jPobbs-Dixon fc Linll2008l : IShowman et al.ll2008b 



Goodmanll2008l ) and we will not attempt to repeat those discussions here but rather briefly 
place this present study in context. In a general sense, the approaches taken to this prob- 
lem can be categorized using two criteria: the approach to radiation and the approach to 
dynamics. Both components are extremely important to include in any model. Energy 
transport, including radiation and advection, determines the distribution of thermal proper- 
ties where as pressur e gradients and gravity regulate thermal currents throughout the planet 
(IBurkert et all 120051 ). A self consistent treatment of both processes is required to determine 
the gas flow pattern and emerging radiative spectra. 

App roaches utilized for radiation to date include relaxation methods (i.e. Newtonian 



heating) (IShowman fc Guif 



Langton fc Laughlinl 120071 . 
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quency dependent radiative transfer (IShowman et al. 
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portion of t he model have includ ed solving the equivalent baro topic equations (ICho et ah 
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Menou fc Rauscher 
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2008 : Rauscher et al. 2008), the shallow water equations jLangton &: Laughlin 



Showman et al. 
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2009 



2007, 



Cooper fc Showman 



20081 : Rauscher Menoi1l20(39[ ) . and Euler's equations ( Burkert et al. 



Ideally we would integrate the full Navier-Stokes equations coupled to a 3D frequency 
dependent radiative transfer model. However, computational resources are limited, and sim- 
ulations of this type take considerable time. Thus, to our knowledge the approach of all 
groups to date has been to choose one approach to dynamics and one to radiation, necessar- 
ily compromising on one or the other. For the results presented here we have chosen to utilize 
a two-frequency, 3D flux-limited model coupled to the full 3D compressible Navier-Stokes 
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equations. Our goal in this study is to understand the details of atmospheric dynamics 
and energy transfer in highly irradiated environments which are not spherically symmet- 
ric. Our new approach to the radiation flux reproduces many of the general features of 
frequency- dependent approaches. Most importantly, it allows for energy deposition at a 
depth corresponding to the optical photosphere. The depth of deposition controls much of 
the subsequent dynamics and a bolometric approach to the energy transfer is all that is 
needed. 

The plan of the paper is as follows; Section ([2]) describes our numeric al approach, concen 



trating on the improvements to our methods from our previous paper (IDobbs-Dixon fc Lin 



20081 ). Primarily, this improvement lies in the treatment of the energy equation in which we 
track the energetic components separately. In addition, we use updated opacities calculated 
for infrared and optical components independently. In Section Q we present results for the 
well known HD209458b. We begin by illustrating the radiative solution, reproducing many 
of the important features of the pressure-temperature profiles from more detailed ID models. 
We then present a number of aspects of the underlying flow and atmospheric structure for 
a range of viscosities. Viscosity plays an important role in both regulating flow velocities 
and providing viscous heating. For a range of viscosities, our simulations exhibit dynamical 
variations that give rise to changing temperature distributions across the photosphere. Sec- 
tion (13. 3p presents several aspects of this exo-weather, in particular changes we expect to 
see from multiple observations. Finally, several questions rega rding the nature of the energy 



flow throughout the atmosphere have recently been posed by iGoodmanI (120081 ). To address 
this, in Section we discuss the role of viscosity in determining the entropy distribution 
throughout the atmosphere, and the relative fluxes of kinetic energy, enthalpy, and radiation 
from day to nightside. We conclude with a discussion in Section ([5]). 



2. Model Description 



We model the planetary atmosphere using a three-dimensional radiative hydrodynam- 
ical model in spherical coordinates {r,(f),6). In addition to the nume rical techniques uti- 
l ized in solving the hydrodynamical portion of the code described in iDobbs-Dixon fc Lin 
( 20081) . we have adde d explicit viscosity terms to our equations in the form described by 



Kley fc Henslerl ( 119871 ). Significant changes to the radiative component o f our atmospheric 
mode l have also been implemented. In addition to updated opacities ( ISharp fc Burrows 
20071 ). we have decoupled the radiative and thermal energy components and introduced a 
multi flux c hannel approach for the rad iative energy component in a manner similar to the 
approach of iHowell fc Greenoughl (120031 ). 



The fully compressible three dimensional Navier-Stokes equation for the motion of the 
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fluid is given by 



- + (u.V)u 



+ g - 2Q, X u - fix (Q, X r) + z/V^u + -V (V ■ u) 

P 3 



where u is the three-dimensional velocity, p and P are the gas density and pressure, and Q 
is the rotation frequency. The gravitational acceleration, g is taken here to be purely radial. 
u = r]/p is the constant kinematic viscosity, proportional to rj, the dynamic viscosity. We 
have neglected the coefficient for the bulk vi scosity. For further de tails on the treatment 
of the viscous terms we direct the reader to iKley fc Henslerl (119871 ). Explicit viscosity is 
supplemented by an additional artificial viscosity for the accurate t reatment of shocks. Th e 
functional form of this term is identical to that in the ZEUS code (jStone fc Normanlll992l ). 
Equation ([1]) is supplemented by the continuity equation, given by 



| + V.(pu) = 0, 



(2) 



The equation for internal energy of the gas can be written as 



g^ + (u.V)e 



-P V ■ u — pen p^i 



Ji 



Era 



(3) 



The quantity proportional to ji represents the photons added to the i-th radiation energy 
density, while the term on the right-hand side proportional to Ep^i represents those absorbed 
by the gas. Kp j is the Planck opacity, defined below. Ignoring non-coherent scattering, these 
terms represent either impinging radiation or interactions with the surrounding gas. In the 
above e = c„/?T is the internal energy density, is the specific heat, and Dy = (SV) v is the 
viscous dissipat i on fun cti on. The components o f the viscous stress tensor S are detailed in 
Kiev fc Henslerl (119871 ) or iMilhalas fc MilhalasI (119841 ). With the exception of the isotropic 
assumption for and (used during the angular integration of Equation (jl]) below) and 
neglecting scattering, this set of equations is completely general and can be expanded to an 
arbitrary number of frequency bins. 

In addition to the internal energy of the gas we must also derive an expression for 
the radiation energy. For this consider a simplified frequency dependent model with two 
frequency ranges. We can write the radiative transfer equation as 



^scat ^scat 



V 



- pi^uh 



(4) 



In the above ly is the specific intensity, is the emissivity, and = k'^^ + k^^"* is the 
total opacity; the sum of the absorption and scattering opacities. The intensity of photons 
scattered into the beam is given by <|)^™* = J I k', k | Jj, (k') dVL\ where (pu ( k',kj is the 
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scattering probability density between the k and k' beams. Integrating Equation (jl]) over 
angle, we define the moments as {cER,i^,Fy) = f (^l^kj Ii,dfl, where Er^i, and F,^ are the 
frequency dependent radiation energy density and flux respectively. 

We now identify the first frequency interval with the incoming stellar radiation. This is 
characterized by the larger T^, and peaks in the optical. The second frequency interval then 
represents the local re-radiated energy, characterized by the local temperature T. This peaks 
in the infrared for the problem of interest here. Addressing the impinging stellar photons first 
(group 1) we assume that the atmospheric gas does not emit in this wave-band, implying 
jl = 0. To find an expression for Ei, the stellar energy absorbed by the gas, we return 
to the full radiative transfer equation. Allowing only absorption, and neglecting temporal 
variations. Equation (jl]) for the intensity of the stellar radiation can be written simply as 

^ = -pK^h, (5) 
ar 

with solution = o exp [— r^,] . We have denoted this frequency interval with the subscript 
-k. The stellar heating term in Equation ([3]) can then be written as 

Er^i = E^ = - [ Ijn = i/..oe"^* = WaTU-^\ (6) 
c J c ' 

where we have explicitly assumed that this radiation is not isotropic but rather flows radially 
from the central star. The factor W accounts for the attenuation of the stellar energy and is 
inversely proportional to the square of the semi-major axis. The optical depth for incoming 
stellar radiation is given by r^,. 

We now return to Equation (jl]) for the local, re-processed radiation (frequency interval 
2), which peaks in the infrared for this problem. Integrating Equation (|4]) over angle and 
frequency, dropping the subscript 2, assuming the gas radiates this energy isotropically and 
is in local thermodynamic equilibrium [i.e. j = Bk), the radiative energy component can 
be expressed as, 

+ V-F = pKp (T, P) [B (T) - cEr] . (7) 



OEr 



dt 

In the above integration we have assumed that ji, and n^, are isotropic. The local radiative 
energy density Er is evolved independently from the thermal component using this equation, 
with the term B (T) = 4(tT^ linking it to the gas. 

Combining Equations ([3]) and iQ, the final thermal energy equation is given by 
de 



,^ + (u.V)e 



-PV ■ u - pnp (T, P) [B (T) - cEr] + pK. (T, P) F.e""* + D,. 



The second term on the left accounts for the advection of thermal energy throughout the 
planet, the first term on the right is the compressional heating term, the term proportional 
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to B (T) — cEji, represents the exchange of energy between the thermal and radiative compo- 
nents through the emission and absorption of the low energy photons, and pn^ (T, P) F^e~'^* 
represents the higher energy stellar photons absorbed by the gas. The stellar flux at the 
planet's dayside surface will be 



[cos (6) cos 



(9) 



where 9 and are latitude and longitude on the planet's surface. 



There are two final ingredients necessary for solving the energy balance in the atmo- 
sphere; a prescription for the opacity and a closure relation linking the flux F back to the 
radiation energy density. Here we utilize the thre e dimensional flux-limited diffusion (FLD) 
approximation of iLevermore fc Pomraningl (jl98ll ). where 



-A- 



-VE 



R- 



(10) 



The Rosseland opacity is given by kr (T, P) , and A is a temporally and spatially variable flux 
limiter providing the closure relationship between flux and radiation energy density. The 
functional form of A is given by 

A^-^i^, (11) 



where 



6 + 3/? + i?2 ' 

NEri 



R= ^ (12) 

FLD has been utilized in a wide range of astrophysical applications and allows for the simul- 
taneous study of optically thick and optically thin gas, correctly reproducing the limiting 
behavior of the radiation at both extremes. In the optically thick limit the radiation energy 
can be expressed as E = aT^, and the flux becomes the standard radiative diffusion equation. 



AacT^ 



-VT. 



3pkr{T,P) 

In the optically thin streaming limit. Equation f lTU]) becomes 

|F| = cEr. 



(13) 



(14) 



Between these lim its, Equation ffTTll approx i mates the full wavelength dependent radiative 
transfer models of ILevermore Sz Pomraning] (jl98ll ). 



There are three frequency averaged opacities appearing in Equations ([7]), (IH]), and ( fTOj) : 
Hp (T, P), (T, P), and kr (T, P). The impor tance of differing absorption opacity and emis- 
sivity was first noted by iHubeny et al.l (120031) but not wide ly appreciated until the discov- 
ery of the stratosphere of HD209458b JKnutson et allboosh . The modeling of HD209458b 
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(iBurrows et al.l 120071 . e.g.) suggests a thermal inversion in the upper atmosphere largely 
due to the differing opacities. Here we include these effects within the framework of fre- 
quency averaged opacities. 3D frequency dependent radiative transfer coupled to the full 3D 
Navier-Stokes equations would be the ideal tool for these calculations; unfortunately, to our 
knowledge no such calculations exist to date. Although our models do not contain the full 
wavelength dependent opacities of one-dimensional radiative models, with the addition of 
multiple average opacities and the radiative energy equation we are able to capture many of 
the essential features that appear in such models. The energetics powering the underlying 
dynamics is well represented by an average approach. Moreover, this approach allows us to 
still solve the full three dimensional Navier Stokes Equations with reasonable computational 
requirements. Our previous atmospheric models JDobbs-Dixon fc Linl[2008l ) not only took 



the radiative and thermal temperatures to be equal but also assumed that all three of these 
opacities were equivalent. Of particular importance, our new approach modifies the radial 
location that the stellar energy is initially deposited and the winds are launched. The un- 
fortunate side-effect is that we cannot produce detailed spectra directly from these models. 
However, it is possible to post-process the pressure-temperature profiles derived here using 
a one- dimensional radiative transfer code, which will be presented elsewhere. 

Planck and Rosseland mean opacities are defined in the usual manner, but the temper- 
ature at which the spectra is evaluated is critical. The local Planck mean for the low-energy 
photon group is given by 

/ (T, P) (T) du 
"^^^'^) = jBAT)du ' ^^'^ 

while the Planck mean for the high-energy photon group is defined as 

>..(r,P)^ J''-gg^(;f'"^ (16) 

The frequency dependent opacity is given by k^, and the subscript ns indicates scattering 
processes are neglected when calculating the wavelength dependent opacities. In principle, 
the impinging radiation can differ from a black-body, but for our purposes we set (T^) = 
By (T^). Finally, the Rosseland mean for the low-energy group is given by 

'"'(r.f)-'--^ aBi,,f ■ (17) 
Here wavelength dependent opacities include the effect of scattering (subscript g) . In the cal- 



culati ons presented here we use the wavelength dependent opacities as described in lSharp fc Burrows 



fl2007f ) 
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Burrows et al.l ( l2008l ) find that matching the spectral data from HD209458b requires 
the addition of a supplemental, unknown uniform opacity source of /tg = O.lcm^/g at pres- 
sures below 0.03 bars. The source of this additional opacity is currently unk nown. Some of 



the l eading contenders are titanium oxide (TiO) and vanadium oxide (VO) ( iHubeny et al. 



20031 ). However, though TiO and VO may indeed be the dominate source of opacity in the 
upper atmosp here at high temper atures (T > 1500K), its abundance remains an outstand- 
ing question. ISpiegel et al.l ( 120091 ) find, in the absence of abnormally large vertical mixing 
coefficients, it is likely that VO s ettles out of the upp er atmosphere and does not play a role 
in producing thermal inversions. IZahnle et al.l ( 120091 ) suggest 5*2 and as another potential 
source of high altitude absorption. Rather then introduce the complex temperature and 
pressure dependence of an unknown amount of Ti O, VO, 5*?, or 5*3, ou r opacities do not 
contain these explicit compounds. Instead we follow iBurrows et al.l (120081 ) and augment the 
Planck opacities in the upper atmosphere with the straightforward Ke parameterization. The 
value of Ke is simply a fit to current observations and is easily modified for studying other 
situations. 



3. Energy Flow Throughout the Planetary Atmosphere 

In this section we present simulations of the full radiative-hydrodynamical solutions of 
Equations ([1]), ([2]), ([7]), and ([8]) for conditions meant to represent the close-in giant planet 
HD209458b. To this end, the orbital and (assumed) rotation periods of the planet are 
P = 3.52days, corresponding to a semi-major axis of 0.0 47 AU. The host s tar is taken 



to have a radius of 1.15Rq and a temperature of 6030K (iMazeh et al.ll2000l ). Although 
observable parameters are chosen from the literature, we have not attempted to tune the 
models to fit observed spectra or light-curves, with the exception of the extra opacity source 
added (see Section (|2])). Our goal here is to explore general properties of atmospheric fiow; 
specific models with tuned parameters will be presented elsewhere. These results can be 
considered representative of fiows as you would see on planets similar to HD209458b. 

We solve the problem in 3D spherical coordinates with a standard resolution of {Nr, N^, Ng 
(60, 160, 64). Due to the decreasing grid size near the pole, and the associated Courant time- 
step constraint, we limit the latitudinal extent of our grid to ±70". These boundaries are 
impenetrable and slip-free, and the fiux is set to zero. The role of these artificial poles on 
the overall fiow is difficult to characterize and is an unfortunate side effect of the chosen 
coordinate system. We have run tests changing \9max\ by 10° and have observed little effect 
on the overall atmospheric structure. However, this test is not conclusive. A more concrete 
discussion of the effect of this boundary requires eliminating this boundary by means of a 
'polar-patch' which we hope to include in future models. 
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All simulations are initialized wi th the static one-dim ensional equilibrium pressure- 
temperature profiles as calculated in iBurrows et al.l ( 120071 ) . As discussed below in Sec- 
tion fl3.ip . we then solve the radiative portion of the energy equations, subject to the non- 
symmetric stellar irradiation until the simulation reaches radiative equilibrium. Dynamical 
results, discussed in Section fl3.2p are initialized with this static radiative-equilibrium model. 
The gas at the inner boundary is assumed to be optically thick [Eji = aT^) and its prop- 
erties spherically symmetric while the boundary is assumed to be impenetrable and slip 
free. We set a spherically symmetric radiative flux and calculate the temperature gradient 
with Equation ( IT3|) . The pressure at the inner boundary is set well within the spherically 
symmetric convective interior, allowing us to both avoid any influence of the choice of inner 
velocity boundary conditions on the upper atmospheric flow and allow for the interaction 
between the convective interior and irradiation induced dynamics. A set of less restrictive 
inner boundary conditions are needed to determine the influence of stellar irradiation on 
convection, cooling, and the quasi hydrostatic contraction rates of hot Jupiters. Such simu- 
lations will be presented elsewhere. A movable outer boundary {not associated with either 
the optical or infrared photospheres) is located at a density of 10~^g/cm^, through which we 
assume a constant outward radiative flux {dFn/dr = 0) through an optically thin gas, given 
by Equation (fT4|) . Given the negligible effect of dynamics at the lowest pressures, the gas 
temperature at the outer boundary is set to its thermal steady-state value with Equation 
( !T8|) . given below. 



3.1. Radiative Solution 

Before including the dynamical portion of the code, we first allow the simulation to 
come to radiative equilibrium. This is accomplished by solving the energy equations given 
by Equations ([7]) and (jH]) with u = 0. The temperature profiles are shown in Figure ([1]). The 
sohd line denotes the gas temperature T and the dashed line denotes the radiative energy 
Eji. Decoupling the energy into thermal and radiative components is a crucial element to 
account for when studying an irradiation induced upper atmospheric temperature inversion, 
and is clearly evident near the stellar photosphere where En remains roughly constant. 
Also shown in Figure ([1]) are the locations of the infrared and optical photospheres and 
the (scaled) direct stellar heating term p/t^ (T, P) -F^e""^* . The extent of the temperature 
inversion is clearly associated with the deeper optical photosphere and the resulting extent 
of the heating term. Ou r radiative pressure-temperature profiles are quite similar to those 



m 



Burrows et al.l ( 120071 ). with the exception of the ad-hoc energy sink they add to the 
dayside at high pressures to mimic the effects of dynamics. This energy sink decreases the 
temperature at depth, producing a secondary inversion in their models. This feature is seen 
self-consistently in our dynamical simulations (Section f l3.2p ). 
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P (bar) 

Fig. 1. — The temperature as a function of pressure for a purely radiative solution (solid 
line) from the sub-stellar point to anti-stellar point. The dashed line shows {Eji/aY^^ and 
illustrates the decoupling of thermal and radiative energies near the photosphere. Solid 
points denote the location of the infrared photosphere, while the square shows the location 
of the optical photosphere on the dayside. Also shown (dot-dash line) is the extent of the 
stellar heating p/t^ (T, P) F^e~^* . 
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To understand the behavior of the temperature at low optical depths consider the steady- 
state behavior of Equation (|H]) for zero velocity. In this case the gas temperature can be 
expressed as 



T 



— + 7 
a 4(T 



1/4 



Local equilibrium between absorption and emission requires that the gas temperature be 
greater then (■^)^^^- In the a bove, 7 is the ratio of Planck opacities evaluated for the stellar 



and local fluxes as defined by iHubeny et al.l (120031 ) as 



7 



(T,P) 

Kp {T,py 



(19) 



The role that the two opacities play in creating a temperature inversion is apparent in 
Equation f|T8|) . the nature of which depends sensitively on the behavior of 7 with depth; 
an inversion requires 7^e~^* to increase outward, while a temperature peak in the outer 
atmosphere (seen in some simulations not presented here) requires a peak in this quantity. 
If the opacity is such that regions exist above the optical photosphere with 7 > 1, the 
gas temperature can easily become super-heated, exceeding the equilibrium temperature 
calculated from Equation (|9]). 



3.2. Mean Dynamical Solution 

In this section we present the full radiative hydrodynamical simulations for HD209458b. 
All dynamical simulations are initialized with the zero velocity radiative equilibrium re- 
sults discussed in Section (13. ip . Fluid velocities are initially artificially damped then slowly 
ramped to their full values to avoid spurious oscillations or shocks associated with the ini- 
tialization. As discussed below, several aspects of the flow and temperature distribution 
have significant time varying components. Results shown in this section are averaged over 
approximately one day, in hopes of representing mean flow properties. Results presented 
here are taken from simulations that have run for a minimum of 100 simulated days. Mean 
temperature and velocity distributions appear to have stabilized by this point. Section (13. 3p 
presents the superimposed time-varying characteristics. 

To explore the role of viscosity we have repeated our analysis for a range of kinematic 
viscosities, shown in Table ([1]). We will refer to runs with different viscosities using the no- 
tation S1-S4. In the absence of more detailed knowledge of the detailed dissipative processes 
{e.g. small scale instabilities or turbulence) we take the kinematic viscosity to be a constant. 
Also listed in Table (JT]) is the characteristic value of ae//,pft = ^/ (cgH) at the photosphere 
for each simulation, where Hp is the pressure scale-height at the photosphere. Although 
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we do not utilize the a-prescription, it is a useful way to normalize our chosen viscosity 
values. Our viscosity is assumed to be ultimately related to the generation of turbulence 
and the cascade of this energy to small length scales. Turbulence may be generated through 
many processes including various shear instabilities (associated for example with overturning 
Kelvin- Helmholtz instabilities), waves acting throughout the atmosphe re (both inertial way es 



in the convective regions and Hough modes in the radiative regions flOgilvie fc Lirul2004l )). 
and shocks generated within the interacting flow. Given the limited resolution achievable 
in dynamical simulations, the cascade of turbulence to scales smaller than our grid size and 
the subsequent dam ping of that energy must necessarily be represented by a sub-grid model. 



Penev et al.l (120081 ) have recently shown that turbulent flows can be well represented by such 
an effective viscosity coefficient, such as in the final two terms in Equation ([T]). Given the 
uncertainty of the sub-grid physics, we have chosen to run models for a range of v. As is 
explained below the values bracket a critical viscosity Vcriti for which we expect viscous effect 
to be (sub)dominate. In addition, there is some amo unt of observationa l motivation for this 



choice from observations of Jupiter's photochemistry. iMoses et al.l (120051 ) show that the eddy 
diffusion coefficient for particles at pressures of 10~'^bars to be 10^ — lO^cm^/s. Presumably 
momentum diffusion will be more effective than species diffusion. 

Figure ([2]) shows the temperature at the infrared photosphere for simulations with a 
range of viscosities. Temperatures range from 1400K to 440K at the photosphere. The 
hottest photospheric point lies slightly east (0 > 0) of the sub-stellar point in all simulations, 
while the coolest points are at higher latitudes (both above and below the equator) on the 
nightside. Along the equator, the location of the coolest point depends on viscosity, lying 
east of the anti-stellar point (0 > vr) for the low viscosity runs S3 and S4. For S2, viscous 
heating associated with the converging flows drives oscillations in the location of the coldest 
point (see Section (13. 3p ). The highest viscosity run (SI) exhibits very little flow on the 
nightside and the coldest point remains very near the anti-stellar point. 



Simulation 


V (cm^/s) 


<^eff,ph 


Hp^ph (km) 


peak V(j, (km/s) 


SI 


10^2 


10-1 


330 


0.8 


S2 


lO^o 


10-3 


360 


4.5 


S3 


10^ 


10-5 


360 


5.3 


SA 


10^ 


10-9 


360 


5.7 



Table 1: Values of kinematic viscosity used for the simulations presented here. For reference 
we also quote an average effective alpha-parameter and pressure scale-height all calculated 
at the photosphere. Peak velocities, shown in the last column, increase with decreasing 
viscosity. 
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Figure ([3]) illustrates the total speed at the infrared photosphere. Here again, the simu- 
lations exhibit a range of behaviors for the varying viscosity, with the highest viscosity (SI) 
keeping velocities below Ikm/s. For the other simulations, as fluid moves from the sub-stellar 
point in the eastward direction [ie to the right in the figures) it is funneled toward the equa- 
tor, while westward moving material is pushed toward the poles. In an inertial frame, planet 
spins from the west toward the east with an angular frequency f2 = 2.1 x 10~^ s~^. Eastward 
and westward winds with 5 km s~^ in the co-rotating frame introduce a local change in the 
azimuthal angular frequency AQ = ±7 x 10~^ s~^ or ±0.3^2. The result is the formation of 
a banded pattern, with jet-like structures. The exact mechanism for the formation of jets 
both in these simulations and in solar system atmospheres remains an outstanding ques- 
tion. The anti-symmetry appears to be associated with the rotation of the planet, perhaps 
rel ated to the asymmet r y of t he Coriolis force (For more detailed plot s of this phenomena 



see 



Dobbs-Dixon fc LinI (120081 )). Non- rotating simulations, presented in iDobbs-Dixon fc Lin 
(120081 ) . show a very symmetric flow pattern despite large wind speeds. Despite the high 
velocity flows, the role of rotation in shaping the flow is not wholly unexpected as the the 
Rossby number, given by the ratio of the advective and rotation terms in Equation ([T]), 
is of order unity throughout most of the atmosphere in our simulations. Given the highly 
unstable nature of the symmetric flow in a non-rotating planet it seems that only a slight 
imposed asymmetry is necessary to form the banded structure seen in Figure ([3]). 

Westward moving material at high latitudes always reaches a stagnation point west of 
the anti-stellar point as it encounters the opposing flow near ~ 140°. At this point the fluid 
subducts under its counter-part, slowing considerably as it crosses down isobars to higher 
pressures and continues around the planet. The behavior of fluid moving in the eastward 
direction at the equator varies with viscosity. For S2, eastward flowing fluid also subducts 
under the opposing flow, much in the same manner as the westward moving fluid. However, 
for S3 and S4 the fluid is able to form a circumplanetary equatorial jet at roughly a single 
pressure level. Peak velocities, increasing with decreasing viscosity to over 5km/s for S4, 
are found just east of the anti-stellar point, but quickly decelerate as it passes back up the 
pressure gradient to the dayside. For reference the sound-speed at the photosphere ranges 
from 1.4km/s on the nightside to 2.9km/s on the dayside, implying the flow velocities can 
reach Mach 3.3 in the equatorial jet. Although S3 and S4 form and maintain a equatorial 
surface circumplanetary jet, there is still significant cooling of the flow after it passes the 
terminator, as can be seen in Figure ([2]). 

To illustrate the subduction of material at the stagnation points, we show the azimuthal 
velocity at the equator {6 = 0) and high latitudes {6 = 35°) as a function of pressure for 
S2 and S3 in Figure (j4]). In both simulations, the wind is launched at low pressures and is 
accelerated as it passes across the terminators (0 = 90° and 270°). These figures also clearly 
illustrate that subduction of fluid occurs for all but the low viscosity simulations at the 
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Fig. 2. — The temperature at the photosphere for simulations with varying viscosi- 
ties. Viscosities for upper-left, upper-right, lower-left, and lower-right are 10^^, 10^°, 10^, 
and lO^cm^/s respectively. The center of each graph is located at the anti-stellar point 
(0, 9) = [it, 0) and the equator (6 = 0) runs horizontally through the center. For lO^^cm^/s, 
the lack of advection keeps the entire nightside cool, while for the other simulations advec- 
tion significantly alters the temperature structure near the terminators (vr/2 and 37r/2) and 
across the nightside. 
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Fig. 3. — The total speed {v = yuf. + u'^ + u^) in cm/s at the photosphere for simulations 
with varying viscosities as in Figure ([2]). Over-plotted are velocities across the photosphere, 
with the arrow length proportional to velocity magnitude. Very little flow is found for high 
viscosity, while supersonic flows are found for the other simulations. The behavior of the 
jets also depends on viscosity, exhibiting both circumplanetary rotation and subduction. 



15 



equator. For example, in the upper left-hand plot of Figure (jl]) the peak velocity of eastward 
moving material (positive v^p) shifts from pressures near 10~^ bars down to pressures near 
10~^ bars near the opposite terminator (0 = 37r/2). The peak velocities of the opposing 
westward fluid remain at low pressures, as the eastward fluid passes under (see the light blue 
profile which changes sign at pressures of approximately 10"'^ bars). The equatorial jet in 

53 and S4 extends over a somewhat larger pressure range than the high-latitude jets or the 
equatorial jet in S2, and also peaks at slightly higher velocities. At high latitudes in both 
S2 and S3 the opposing streams meet closer to the terminator east of the substellar point 
(0 = 90°). Although some degree of subduction is seen in the yellow profiles, most of the 
fluid is pushed toward the equator, where significant down-welling occurs. The flow toward 
the equator can be seen near </> = 90° in Figure ([3]). 

The transition in the overall flow structure (subduction vrs circumplanetary jet) between 
S2 and S3 can be related to the critical value of Ucrit for which viscous effects become 
comparable with advection. To quantify this we calculate a modified Reynolds number, 
where the length scales used for calculating the viscous and advective terms are different. 
We express the crossing timescale as tx = Rp/uavg and the viscous timescale as t^, = Hp/u. 
We use the vertical scale-height Hp because we find in our simulations that the radial transfer 
of momentum (oc ^^) is overwhelmingly the primary viscous sink for the jets as they traverse 
the planet. The ratio of r^, to gives our modified Reynolds number. Re. Averaging over 
the photosphere we find i/p ~ 3 x 10'' cm, Uavg ~ Ikm/s, and i?p ~ 9 x lO^cm. Setting 
Re = 1, we find that Ucru ~ 10^°cm^/s is the critical viscosity, above which viscous forces 
begin to dominate. For u > Ucru viscous forces try and maintain the symmetry in the flows 
from both east and west. The flow patterns for S2 seen in Figure (jl]) represents the transition 
between advection dominated (simulations S3 and S4) and viscosity dominated (simulation 
SI) flows. 

In Figure ([5]) we show the pressure-temperature profiles at the equator for all four 
simulations. The inversion at low pressures (< 10~^ bars) and interior adiabat are both still 
clearly evident as seen in the radiative profile shown in Figure ([1]). Temperatures extend 
up to 2600K near the inner radial boundary, where the planet is spherically symmetric and 
fully convecting. Dayside temperature profiles (shown in red) are fairly consistent across 
all simulations S2-S4. The upper regions of SI remain slightly hotter, with temperatures 
closer to the original radiative profile of Figure ([1]). The larger velocities calculated in S2- 

54 are more effective at advecting away energy, cooling the dayside. However, the relative 
similarity between the dayside photospheric temperature of S2-S4 suggests that the profile of 
infrared and optical opacities (Equation (fT9|) ) is the dominant factor in determining dayside 
temperatures at pressures lower than 10~^ bars. 

The major effect of the fluid flow on the subsolar profile is at higher pressures, below both 
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Fig. 4. — The azimuthal velocity [cm/ s) as a function of pressure at the equator for S2 (top 
row) and S3 (bottom row). The left column is a slice at the equator [9 = 0°) and the right 
column is a slice at mid- latitudes [6 = 35°). Profiles proceed from red at the sub-stellar 
longitude all the way around the planet in 18° steps in longitude. The key in the upper-left 
panel shows the longitudinal location of a given profile looking down the rotation axis, with 
the arrows indicating the direction of incident stellar flux. Solid dots denote the location of 
the infrared photosphere, while the square shows the location of the optical photosphere at 
the substellar point. 
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the optical and IR photospheres. The returning flow significantly cools the atmosphere at 
pressures around 1 — 10~^ bars creating a secondary, dynamically induced inversion at depth. 
The need for two inversions in the pressure/temperature profile, one caused by the different 
absorption and emission op acities and the other caused by the dynamics, was introduced in 



a parameterized manner by iBurrows et al.l (120071 ) in order to simultaneously match spectral 



constraints at multiple wavelengths. It arises self-consistently here. Comparing Figures (jlj) 
and this cooler region is associated with the flow returning from the nightside after it 
has circumnavigated the planet and cooled while on the nightside. The final equilibrium 
temperature at this location will be determined by a balance between the incident energy 
diffusing downward from the optical photosphere and the advection bringing in cooler fluid 
from the nightside. For this reason, this secondary inversion is slightly warmer and narrower 
(in pressure) in S2, than in S3 or S4. As seen above, the higher viscosity (S2) limits flow 
velocities, giving the downward diffusing radiation more time to heat the fluid. 

The differences between S2-S4 are most prominent at regions near the terminator east 
of the anti-stellar point (see the darker blue profile in Figure ([5]). The deep, fast jet of S4 
cools the terminator east of the anti-stellar point (</> = 37r/2) to a much greater extent then 
in any of the other simulations. This cooling produces a very distinct temperature inversion 
at this latitude. Such inversion features should be easily distinguishable in a spectra taken 
after primary transit but before secondary eclipse, when this terminator faces toward earth. 
In contrast, the terminator west of the anti-stellar point exhibits no inversion in any of the 
simulations. 

The inclusion of viscosity in the governing equation provides another mechanism for con- 
verting kinetic energy back into thermal energy, regulating the velocity in the process. The 
left-hand plot of Figure illustrates the strength of the viscous term [i^V^u + |V (V ■ u)] ^ 
in the longitudinal momentum equation for S2. This plot concentrates on the convergence 
region on the nightside, and clearly illustrates that viscosity primarily acts in this shear zone. 
Comparison to Figure @ shows that the two shearing streams have relative velocities of up 
to 8km/s. The heads of the jets experience the greatest deceleration, reaching 1250cm/s^. 
The importance of viscous dissipation on the overall energy distribution depends on its con- 
tribution to Equation (|8]). The other component making a major contribution to the thermal 
energy is the compressional work, PV ■ u. Though for a majority of the atmosphere PV ■ u 
dominates, Di, becomes more important in the convergence region for the lower viscosity 
simulations. In contrast, for the highest viscosity simulation SI, viscous heating is most 
relevant near the terminators where the velocities are highest. The right-hand plot shows 
the ratio of the viscous and compression terms on the nightside near the anti-stellar point. 
D^, exceeds PV ■ u by several hundred times at the jet heads, which extend over a wide 
range in longitude. The details of this convergence region are paramount in determining 
the overall dynamical structure throughout the atmosphere, as it is here where flow either 
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Fig. 5. — Temperature (solid lines) and radiation energy density (Eji/a)^^^ (dashed lines) 
profiles as a function of pressure at the equator for S1-S4. Colors and symbols are as in 
Figure (jH). 
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subducts under its counterpart or instigates circumplanetary jets. Neglecting the viscous 
contribution, especially in this region, may lead to erroneous results. Across all simulations, 
we find that the average ratio of to PV ■ u near the photosphere roughly scales with the 
value of the kinematic viscosity. 

The relation between nightside temperature, velocity, and viscosity is important to 
understand. Many authors have noted that the nightside temperature can be understood by 
equating the crossing timescale with the radiative cooling timescale; larger velocities carry 
more energy to the nightside, decreasing the temperature differential, and thus the azimuthal 
acceleration. It is easy to see an equilibrium will quickly develop. The introduction of 
viscosity complicates the issue, controlling the velocity and heating the gas. Substantial 
viscous heating requires both shear and large v. In SI, the high viscosity prevents fiow to 
the nightside and thus viscous heating is confined to the dayside. Near the terminators 
exceeds PV ■ u, but the energy is still dominated by the incident stellar heating, and the 
viscous contribution is insignificant. In contrast, there is significant shear on the nightside 
of S4, but, due to the low z/, D^, is an order of magnitude smaller. S2 seems to be in a sweet 
spot. The viscosity is low enough to allow for shear fiows to develop on the nightside, but 
high enough to still generate significant heating there. 



3.3. Time Dependent Behavior: Exo- Weather 

Although we have presented average quantities in the Section (13. 2p . the atmosphere can 
have a dynamic component super-imposed on the stationary background state. As described 
above, the background state consists of a quasi- stationary fiow pattern with fiows traveling 
in both eastward and westward directions. As a result, sustained supersonic winds traveling 
in opposite directions converge on the nightside in the upper atmosphere. Shocks, shearing, 
and instabilities result, and the location of the convergence point oscillates over a range 
of latitude and longitudes. The shock heating {Dy) and compressional heating (PV ■ u) 
associated with this convergence causes the temperature distribution on the nightside to 
fiuctuate on the timescale of approximately half a day. Figure ([7j) shows the temperature 
distribution at the photosphere over a total of 1.5 days. As can be seen, fiow variations 
drive temperature changes on the nightside, with the largest deviations occurring near the 
terminators. 

Changes in the velocity structure are driven primarily by the interaction between high 
velocity shearing fiows. Shearing occurs both between streams at constant pressure and also 
between different pressure layers. Shearing at constant pressure can be seen on the nightside 
in Figure ^ at ±15° latitude for S2-S4 where eastward and westward material pass by each- 
other and interact. Shearing between flows at different pressure levels can be seen in Figure 
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(jl]), where the subducting material travels under counter- flowing material. This can be seen 
at ^ = 35° for S2 and S3, and additionally at the equator for S2. We flnd that variations in 
the temperature at the photosphere are principally driven by the subduction at the equator, 
as eastward moving material sinks below westward moving material at cj) ~ 37r/2. Signiflcant 
potential vorticity is generated in this subduction region, with vortices forming between the 
two shearing flows. As the convergence point moves across the surface it forces these vortices 
to interact, stirring temperatures further, both at the photosphere and at depth. Over the 
course of the entire simulation the location of the coldest point in the calculated phase-curve 
is seen to oscillate more than 20° in longitude, occasionally shifting to the west of the anti- 
stellar point. The largest temperature variations (up to 15%) are seen along the terminators, 
where changes to flow patterns have the largest impact on advective efficiencies. 



4. Enthalpy, Radiation, and Kinetic Energy Fluxes 

One important, overarching questions for groups that model the atmospheric dynamics 
of hot-Jupiters is calculating how much energy is transferred to the nightside, and how 
changi ng conditions ch ange the resulting dynamics and energy distribution. In a recent 



article iGoodmanI ( l2008l ) points out that without viscosity the atmosphere will adjust to a 
different solution in order to achieve a steady state with zero Carnot efficiency. This is 
fundamentally an argument regarding the momentum equation, not the energy equation. 
To illustrate this, consider the momentum equation ([1]) in steady state, neglecting rotation. 
Taking the dot product with velocity and expressing the advective term as (u ■ V) u = 
1/2 V |u|^ + (V X u) X u, we can write 

u ■ V Q |u|^ + w + = u ■ TVS + u • z/V^u + u ■ (V ■ u) . (20) 

The term in parenthesis is the Bernoulli constant {Eb) of the flow, S is the entropy per unit 
mass, and w = e/ p + Pj p is the enthalpy per unit mass. Bernoulli's constant is the total of 
the kinetic, thermal, and potential energies. Next consider the steady state form of the two 
energy Equations ([7]) and ([8]), which together become 

(u ■ V)e = -PV ■ u - V ■ F + + D^. (21) 

Si, = pKi, (T, P) F^e""^* is the direct stellar energy absorbed by the gas. Utilizing the flrst-law 
of thermodynamics, this can be expressed as 

u ■ TWS = p^^[D^-W ■¥ + Si] , (22) 

where the PV-u term cancels out, as this does not contribute to changes in entropy. Plugging 
this into Equation (120|) we see that the gradient of Eb along streamlines is 



u ■ VEb = p^M^i^ - V ■ F + + u ■ z/V^u + u ■ ^ V (V ■ u) . (23) 
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However, the role of viscosity here is to convert kinetic energy to thermal energy. Therefore, 
the viscous contribution from the momentum equation must be equal and opposite to the 
viscous contribution from the energy equation {p~^Dj, = — u ■ z/V^u — u ■ |V (V ■ u)) leaving 
us with 

u ■ VEb = p-' [5. - V ■ F] . (24) 

Changes in the gradient of Eb along a streamline are solely due to radiation effects, even when 
viscosity is included. Since the integral of u ■ VEb around a closed streamline must equal 
zero in steady-state, so too must the radiation terms. Therefore, the important difference 
when including viscosity is the behavior of the streamlines. Including viscosity will alter 
the velocity, changing the streamlines and thus also the overall thermal structure of the 
atmosphere. Figures ([3]) and (jlj) show that viscosity can play an important role in altering the 
velocity (and thus thermal) structures for the planetary atmosphere. Ultimately, Equation 
(!24|) tells us that the direct effects of viscosity are important for solving the momentum 
equation, not the energy equation, further justifying our choice to solve the full 3D Navier- 
Stokes equations. 



As pointed out by iGoodmanI (120081 ). another way to understand the importance of vis- 
cosity is to calculate the relative importance of the longitudinal energy fluxes carried as 
kinetic energy vs that carried as enthalpy. The left-hand panel of Figure ([8]) compares the 
azimuthal flux of enthalpy {e/ p + P/ p) to the azimuthal flux of kinetic energy (|u|^ /2). En- 
ergy carried in the form of enthalpy dominates that carried in kinetic energy form throughout 
most of the interior. However, at higher altitudes the two contributions become comparable. 
In particular, in regions of high velocities on the nightside, where we have already seen that 
shocks play a role in determining the dynamical structure, the flux of kinetic energy exceeds 
that of enthalpy. In these regions energy is primarily carried in kinetic form and energy 
loss through shocks is important. Our approach of integrating the Navier-Stokes equations 
allows us to follow the transfer energy from kinetic to internal energy across sh ocks even if 



the shocks themselves are not fully resolved. As emphasized by lGoodmanl (120081 ). this in not 
the case for more simplified treatments of the dynamics. 

One final advantage of the method we are using here is that radiation is allowed to 
flow in all three dimensions. Incident radiation is modeled here as being purely radial, but 
re-radiated photons can also carry energy both to the nightside and pole-ward. Although 
the radial component is by far the most important, at high altitudes with pressures less 
than ~ 0.01 bars, the azimuthal flux of radiation can become comparable to the radial 
radiative flux. The ratio of F^ to F^ to is shown in the right-hand plot of Figure ([H]). 
Although never exceeding the radial radiative flux, the azimuthal radiative flux is of the 
same order throughout the dayside and extending to the terminators. The effects of the 
azimuthal radiative transfer are most prominent near the terminators, where they play a 
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role in shaping the final pressure-temperature profiles. 



5. Conclusion 

We have performed a series of simulations for the atmosphere of the highly irradiated 
planet HD209458b. Our approach involves solving the full 3D Navier Stokes equation cou- 
pled to both the thermal and radiative energy equations. We have separated our frequency 
averaged opacities into those relevant for incoming stellar photons and those for re-radiated 
photons. This approach allows us to reproduce the observed temperature inversion in the 
upper atmosphere and self- consistently calculate the location that stellar energy is deposited. 
This deposited energy drives fiow in both the eastward and westward directions, with east- 
ward fiow concentrated near the equator and westward fiow pushed to higher latitudes above 
and below the equator. These fiows then converge on the nightside and exhibit range of be- 
havior for various explicit turbulent viscosities. 

The viscosity in our simulations spans a wide range, running from v = 10^^ to lO^cm^/s, 
corresponding to effective ae//,pfc-parameters ranging from approximately 10~^ to 10~^. As 
expected, results vary significantly with viscosity, especially for higher viscosities. For our 
highest viscosities, atmospheric velocities remained subsonic, advection of energy is sup- 
pressed, and the day-night temperature contrast remains high. As viscosity is decreased 
fiow velocities became super-sonic, reaching speeds of over 5km/s in the lowest viscosity 
simulation. In addition to the speed, the behavior of the jets significantly changes with vis- 
cosity. For the simulations with low viscosity {aeff,ph = 10~^ and 10~^) a jet develops that 
circumnavigates the planet at a single pressure near the equator. However for intermediate 
viscosity (ae//,p/i ~ 10"'^), and at high/low latitudes in all simulations, circumplanetary fiow 
of a jet can only be achieved accompanied by subduction. Subduction occurs as the fiuid 
cools and sinks radially, crossing isobars. 

Viscosity also plays a role in determining the variability of the atmospheric fiows and 
energy distribution. This variability manifests itself in two ways: temporally varying ve- 
locity structures that alter advection efficiency, and through viscous heating. Seen most 
prominently in the simulation with u = lO^^cm'^/s (S2), shearing between converging fiows 
on the nightside in both the horizontal and radial directions drives oscillations in the fiow 
patterns. These oscillations alter the location of the convergence point and shape of the 
streamlines thus changing the rate of energy advection coming from the dayside, primarily 
affecting the temperature near the terminators on the level of 15%. Changing streamlines 
also shift the region of maximum shear, in turn altering the location of hot-spots on the 
nightside induced by viscous heating. The location of the coldest point in the calculated 
phase curve can shift by more than 20° in longitude, occasionally even appearing west of 
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the anti-stellar point. Such time-dependent behavior is somewhat suppressed in our other 
simulations (SI, S3, and S4). Viscously driven variations in the flow require that the viscous 
term (z/V^u) in the momentum equation be of order the advective term ((u ■ V) u). This 
only occurs for sufficiently high shear (present for low-z/ simulations) and also large enough 
kinematic viscosity [u). The u = lO^^cm^/s simulation exists in somewhat of a sweet-spot, 
allowing high velocity flows and still maintaining significant viscous heating on the nightside. 
Such variation, if observed, may provide a clue to the nature of the turbulent viscosity in 
close-in giant planet atmospheres. 

Fluctuations in velocity and temperature are largely confined to the nightside. Periodic 
dynamical variations on the dayside are not able to alter the basic structure at an observ- 
able level due to the dominant incident stellar flux. This is a somewhat different result then 



the si mulations of ICho et al.l ( l2003l . l2008l ) ; iRauscher et al.l ( 120081 ) and iMenou fc Rauscher 
( 120081 ). where they demonstrate significant variations throughout the entire phase-curve for 
a range of selected parameters. The variations in their models are largely due to cold cir- 
cumpolar vortices, though they do observe smaller amplitude variations near the equator. 
We do not see such cold circumpolar vortices, though their development may be impeded 
due to our longitudinal boundary conditions. Variations driven by large scale vortices are 
considerably different in nature then those presented here, which arise from the shearing 
between converging flows and are confined primarily to the nightside. We find, apart from a 
modest downwind advection of the hot-spot, dynamics does very little to alter the temper- 
ature profiles very high in the atmosphere on the dayside. Simulations with high velocity 
winds cool the upper atmosphere slightly, but the final temperatures are fairly similar to the 
radiative simulations. The pressure-temperature profile at pressures less than 10~^ bars is 
thus prim arily dependent on o ur choice of opacity and its behavior with depth as originally 
shown by iHubeny et al.l (120031 ) . With a lack of better knowledge concerning the abundance 
of high opacity sources such as TiO or VO, we have chosen to augment the opacity with an 
ad- hoc extra s ource at high pressures (IBurrows et al.ll2008l ). For comparison, the models of 
HD209458b by lShowman et al.l (120091 ) include TiO and VO in the calculation of opacity, and 
as a result of this larger (T, P) achieve much higher temperatures at high altitude. Given 
that this extra opacity is the dominant factor in determining the temperature, matching 
future observations for a wide range of temperatures requires only adjusting the magni- 
tude of this opacity . Tera poral variation in this opacity would drive variations in dynamics, 
though lAgol et al.l ( 120081 ) studied a series of secondary eclipse of HD189733b and find that 
a variation in F^ay/F^ of less then 10%, consistent with a stable dayside temperature. Tran- 
sit spectroscopy appears to be the most promising observational technique for constraining 
exo-weather. 

The implications of dynamical behavior across the nightside is clear; it changes the 
temperature and location of the coldest spot, and it changes the temperature structure 
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near the terminators. The concentration of the fluctuations near the equator serves to 
enhance the observabihty of such variations. Average nightside temperatures, important for 
constraining interior and evolutionary models, can only be obtained by averaging a number 
of such observa tions. This re s ult al so suggests a possible interpretation for the puzzling 



results seen by iKnutson et al.l (120071 ). in which they observe a local minimum in the flux 
before the secondary eclipse for HD 189733b. Though we have no simulated this planet, in 
light of the dynamical behavior seen in our simulations the increase in flux before secondary 
eclipse may be attributed to transient exo-weather, which shifts the shape of the planets 
phase-curve over time. 
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Fig. 6. — The azimuthal component of the viscous force at the photosphere (left), acting 
on U0. The right panel shows the ratio of viscous and compression contributions, plotted as 
loQio (l-Dy/PV ■ u|). Both graphs are contoured with the total velocity at the photosphere 
and are for S2, with u = 10^°cm^/s. We have zoomed into the region on the nightside 
centered on the anti-stellar point with the equator running horizontally through ^ = 0, to 
better illustrate this crucial region. 
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Fig. 7. — The time dependence of S2 {aeff,ph = 10~^). Snapshots are shown half a day apart 
and illustrate temperature (color scale) and total velocity (contours). Radial and horizontal 
shearing drives variations in the flow, changing the advection of energy across the terminators 
and influencing the location of viscous dissipation. 
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Fig. 8. — The ratio of enthalpy and kinetic fluxes throughout the photosphere (left panel). 
In regions with Fk > Fe, it is important to include viscous effects. The right panel show the 
ratio of the radial to azimuthal radiative energy flux. Azimuthal radiative transport becomes 
particularly important near the terminators. Both plots are shown as the logarithm of the 
absolute value of the flux ratios and are for S2. 
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